Model for the hydration of non-polar compounds and polymers 
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We introduce an exactly solvable statistical-mechanical model of the hydration of non-polar 
compounds, based on grouping water molecules in clusters where hydrogen bonds and isotropic 
interactions occur; interactions between clusters are neglected. Analytical results show that an 
effective strengthening of hydrogen bonds in the presence of the solute, together with a geometric 
reorganization of water molecules, are enough to yield hydrophobic behavior. We extend our model 
to describe a non-polar homopolymer in aqueous solution, obtaining a clear evidence of both "cold" 
and "warm" swelling transitions. This suggests that our model could be relevant to describe some 
features of protein folding. 
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I. INTRODUCTION 

The physical properties of liquid water, despite the 
extensive studies devoted to their investigation, are not 
yet completely understood. The behavior of water when 
nonpolar solutes are present is still a matter of debate, 
even if simplified models [0|2j as well as numerical simu- 
lations (see Ref. and references quoted therein) have 
shed some light on it. The hydration of non-polar com- 
pounds has a free energy cost, resulting in a net "attrac- 
tion" among solute molecules which tend to aggregate 
in order to minimize exposition to water ("hydropho- 
bic effect"). More precisely, upon hydration of a sin- 
gle molecule from the gaseous phase, the free energy 
and specific heat changes AF and AC are positive over 
the entire temperature range accessible to experiments; 
the internal energy end entropy differences AE and AS 
are negative for low temperatures, then become positive 
when the temperature is raised. The two temperatures 
at which AE and AS are zero are commonly referred 
to as Th and Ts, respectively. The origin of the lower- 
ing of the entropy for T < Ts has been related to the 
ordering of water molecules around the solute to pre- 
serve bonds [^,^ as well as to the opening of a "cav- 
ity" in water to settle the solute, without ordering of 
water molecules Q; but the question is still controver- 
sial. Hydrophobicity plays a fundamental role in many 
physical processes, among which the long-standing prob- 
lem of protein folding, since the hydrophobic effect is 
believed to be one of the most relevant "forces" guiding 
the protein into the folded native state, i.e., the unique 
spatial structure in which the protein is biologically ac- 
tive j^. At present, the whole folding process can be 
studied only using coarse-grained models (all-atom sim- 
ulations of proteins and solvent molecules are feasible 
only on timescales several order of magnitude smaller 



than the folding time). In many models water is not 
explicitly considered and solvent-induced aggregation is 
described through temperature-independent effective at- 
traction among nonpolar monomers. An important fea- 
ture of real (globular) proteins is that the native state be- 
comes unstable not only at high temperature (warm un- 
folding), but also at low temperatures (cold unfolding) [Bj; 
cold swelling has also been detected in some homopoly- 
mers |9|,[l0t]. A major drawback of the above simplified 
models is that they cannot account for cold unfolding, 
because there the native (folded) state is identified with 
the model ground state, so that the protein is folded at 
low temperatures down to T = 0. 

In this paper we propose a "minimal" model capable of 
reproducing the thermodynamic hallmarks of hydropho- 
bicity in single molecule's hydration, within the frame- 
work of the "water ordering" hypothesis; moreover, we 
keep it as simple as to be analytically integrable, in or- 
der to apply it to polymer hydration. For its simplicity, 
the model is not intended to reproduce all the details of 
hydrophobic behavior: rather it gives a correct overall 
description of the phenomenology and allows to under- 
stand what are the key ingredients that produce it. 

In the case of a nonpolar homopolymer we recover both 
"cold" and "warm" swelling transitions, thus strengthen- 
ing the idea, already put forward in Refs. that an 
explicit, though simplified, description of water can also 
provide a framework for a unified treatment of both tran- 
sitions in proteins. 



II. THE MODEL 

We describe water molecules as two-dimensional ob- 
jects with three hydrogen-bonding arms, namely, the pla- 
nar projection of the tetrahedral coordination of water 
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molecules, as in Refs. [pl3|Jq]. Simulations show that such 
a crude model reproduces the main features of water . 
The molecules in contact with a solute will be referred 
to as the hydration shell. Hydrophobicity is a collective 
phenomenon, yet experiments |^ show additivity in the 
nonpolar area exposed to the solvent. Hence, there ex- 
ists a cluster scale above which correlations may be ne- 
glected. We consider the hydration shell to be made up 
of z such clusters, each of m interacting molecules, and 
we disregard interactions between clusters. We assume 
that the ground state of a cluster is characterized by a 
completely formed hydrogen-bond network, both in the 
bulk and in the hydration shell case. Thus each molecule 
in the ground state has on average 3/2 hydrogen bonds. 
Hydrogen-bond energies are very sensitive to the geom- 
etry of the bond, and the latter will be affected by the 
presence of a solute, so that we assume that bond ener- 
gies may differ in the "bulk" and "shell" cases. We model 
this by simply assuming that there is an energy difference 
(per molecule) K between bulk and shell ground states, 
given by 

K^^{hi^hs) + J , (1) 

where /if, and kg are the hydrogen-bond energies for bulk 
and shell molecules, respectively. The fact that /if, 7^ 
is not, in principle, the unique possible cause of an en- 
ergy shift in the ground states: not only a solute perturbs 
the water hydrogen bonds, but also has a different direct 
interaction with water. To take into account this kind 
of contributions to the ground-state energy difference, 
assuming that they can be reasonably averaged on the 
molecular scale, we introduce the constant J in Eq. (|l|) 
above. 



A. Partition function of the hydration shell 

Our goal is to model the thermodynamics of hydration, 
hence we need to evaluate the partition function for the 
cluster in both cases: 

/■oo 

/ de5.(e)e-'^^-'"("), (2) 

where • = 6, s (bulk or shell, respectively). Here, 

= £ + Xm(5.^, , (3) 

and g,{e) is the density of states of the cluster at the 
energy e above its ground state: it encodes all the im- 
portant features determining the system thermodynam- 
ics. How can we estimate (7, (e)? Given an energy e to 
the cluster we expect to have a certain number of broken 
bonds and assume equipartition of the energy on the var- 
ious degrees of freedom of the cluster. Moreover, bond 
breaking allows molecules to access a wider configura- 
tion space: while the ground state is unique, different 



geometrical arrangements of the molecules can have the 
same number of bonds. To describe it simply enough, we 
take each bond to be represented by a harmonic poten- 
tial, and to be broken independently of the others. Thus, 
at the cluster energy e there will be a probability 

p,{\,s) ^ (1-Pi.(£))^'"-' (4) 

to break A bonds. Here pi,(e) is the probability that 
a bond acquires an energy larger than /i, and breaks, 
when e is equipartitioned on D degrees of freedom: 
pi,{£) — e^^'l^ and T = lejD. For the sake of simphc- 
ity we assume D — 2m f (/ are the degrees of freedom of 
one molecule: / = 3 in two dimensions) as if molecules 
were at least singly bonded at all energies. In the same 
spirit we take the density of states oj{X, e) of the system of 
harmonic oscillators, resulting when A bonds are broken, 
to be independent of A: uj{X,e) ~ oj{e) — Ce~^^, where 
C is a constant. To account for the geometric degeneracy 
we assume just two orientations for each molecule, (E and 
I; see Fig. and give an estimate of 7, (^), the number 
of accessible geometric arrangements of the cluster when 
A = 3771/2 — 1/ bonds are broken. Then, our ansatz for 
5, (e) will be the sum, over A, of the probability p,(A,e) 
of breaking A bonds, times the density of states uj{e) of 
each geometric configuration, times the number "fm{v) of 
configurations at A broken bonds: 

In the estimation of ^m{v) the bookkeeping is easier if 
we attribute each bond to just one of the two bonding 
molecules. Each cluster can have no more than 3m/2 
bonds: to ensure this, we take m even and fictitously 
split the cluster in two equal sets, a and b. Group a- 
molecules will form at most one bond: a shell one when 
in state-E and an external one when in state-I. Group 
^-molecules form at most two bonds: the external and 
a shell one in state-E, an external and the internal one 
in state-I. Internal bonds will be possible only for (&,I) 
molecules, if they find a bonding partner. The number 
of internal bonding sites, q,, will distinguish "shell" from 
"bulk" cases: in the former Qs = 0, while in the latter 

< Qb < m. In fact, due to the geometry of the internal 
water molecules, possibly not all the m positions allow 
internal bonds, even in the bulk case. State-E molecules 
can form one shell bond, while state-I cannot: hence, the 
probability of a configuration with s state-I molecules 
separating (m — s) state-E ones, so that k shell bonds are 
forbidden and x = ^ti — s — fc can be formed, is 

The probability that one of these configurations also has 

1 internal bonds depends, first, on the probability of fish- 
ing out j (6,I)-molecules among the s in state-I, given 
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the total number of molecules m and of 6-molecules m/2: 
7rh(j; s, to/2, to), where is a hypergeometric probabil- 
ity, i.e., 

-(—)^ (!)(::!)[(:)]"■■ (^) 

Then, one must consider the probability of placing i of 
the j molecules with internal arms, in the gr, good po- 
sitions for bond formations, on a total of to possibili- 
ties: TThii] j , q» , m) . The product of all the above prob- 
abilities gives the fraction of conformations that has s 
state-I molecules and is able to form x shell bonds, i 
internal ones and to/2 -t- s — j external ones (assuming 
that all the external arms form bonds). If we now let 
the bonds be also broken, we must notice that any pat- 
tern geometrically allowing more than v bonds will also 
contribute to 7»(i'): we just need to consider the ex- 
cess bonds as broken - we choose them among the exter- 
nal, internal and shell ones with a binomial probability 
TTb {n, N) = TTb {n, TV, 1/2), where 

Mn,N,p)=(^^^p'\l-p)^''--K (8) 

Hence, upon summing over all the geometric arrange- 
ments that can contribute to a pattern with h' bonds, we 
finally get 

1/^=0 i/e=0s=0fe=0j=0i=0 

where Vi — v — Ug ^ v^, fcmax = min(s, m — s — Ug) and 

(jn \ 
r,s,—,mj 7rh(i;j, q„m) 

XTTb {Vi, i) TTb [Ve, — + S - j j TT], {vs,m - s ~ k) . (10) 

The partition function can now be evaluated from Eq. 
(|) and Eqs. (||,|). We get: 

Z,^'" = ^.e-^^"*''= (11) 

where A, = C (\") 7. (3to/2 - A) /. (A) , C is a 

constant and 



^(A) = 




with Kd{t,) the Bessel-X function, A' = 3to/2 — A, and 
a, = [213)-^ Dh, {j -I- A) , T. = ^2(3h,D (j + A) . (13) 



B. Statistical mechanics of non-polar hydration 

We now consider the transfer of a nonpolar molecule 
from its gas phase to water (i.e., hydrophobic hydration). 
We describe it by taking z clusters, of to water molecules, 
in the bulk case and substituting the internal water with 
the solute. The free energy change A_F is given by 

/3Ai^ = -zlogx, (14) 

with X = AsA^"^ exp{—/3Km). Similar expressions hold 
for energy, entropy and specific heat changes. The tem- 
perature dependence of these functions, reported in Fig. 

shows the hallmarks of hydrophobic behavior: we find 
a maximum in the free energy cost and minima in both 
AE and TAS; as T grows, first AE and then TAS cross 
the zero - at Th and T5, respectively - and eventually be- 
come positive. Moreover, we correctly find a pronounced 
and positive peak in the specific heat difference AC. This 
cannot be found within the information theory approxi- 
mation to the cavity model: see Ref. [|l5|. It turns out 
that both a ratio hs/hb > 1 and the difference between 
7b (i') and Jsiv) concur to produce the correct thermo- 
dynamic behavior: with = hh the difference between 
7b (z/) and 7s (z^) is not sufficient to produce the minimum 
in AE, while 76 (z^) ~ ls{v) (i-C, same chance of internal 
bonds, qb — qs) yields too small a peak in AC. 

Then, according to our model, hydrophobicity ap- 
pears as related to the average strengthening of hydrogen 
bonds, together with geometric pattern changes, upon 
the introduction of a nonpolar solute. 

In particular, having tighter shell bonds than bulk ones 
{hs > hb) appears to be essential, in our model, to get 
a negative minimum in AE. Indeed, even starting from 
equal energy ground states in the bulk and shell cases 
(i.e., K = 0), it hs > hb it will be easier to break bonds 
in bulk than in shell, according to Eq. (^. This yields 
that at "intermediate" temperatures shell water will have 
lower energy than bulk, i.e. AE < 0. 

The predictions on the trends of thermodynamic func- 
tions are in qualitative agreement with experimental re- 
sults for noble gases (see the results by Crovetto and 
coworkers [|l6| as reported in ||]) and for nonpolar sur- 
faces in proteins ||^ . They are also in agreement with the 
Monte Carlo results for the two-dimensional "MB" wa- 
ter model in Ref. ||^, where an average strengthening of 
shell hydrogen bonds is detected, too, as well as negative 
values of AC at low temperatures. 

It should be mentioned, though, that negative AE and 
AS*, and positive AF can be obtained for the MB model 
within the information theory approximation applied to 
the cavity opening in pure water Ref. |l^, even if, as 
already stated, this approach cannot recover the correct 
trends for the specific heat. Since cavity opening disre- 
gards bond differences in bulk and shell, this fact would 
imply that hydrogen bond differences are irrelevant to 
the subject. 
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It is likely that both cavity formation and hydrogen- 
bond differences contribute to the hydrophobic hydra- 
tion, and an ultimate answer will come from simula- 
tions in three dimensions when they will become accurate 
enough as to predict specific heats. As far as the present 
model is concerned, we consider the good overall agree- 
ment with experiments and two-dimensional simulations 
as satisfactory to justify the use of this model in polymer 
hydration studies. 



C. Cold and warm unfolding of a non-polar 
homopolymer 

Let us now turn to the study of a nonpolar homopoly- 
mer in solution, taking into account just the behavior of 
water clusters in the vicinity of a monomer, and disre- 
garding interactions between monomers and between wa- 
ter clusters. Our goal is, in fact, to understand the effect 
of the hydrophobic effect alone on polymer behavior. We 
model a polymer as a A^-step self-avoiding walk (SAW) 
on a two dimensional lattice with coordination z. On 
each lattice site there can be either a monomer or z clus- 
ters of m water molecules, so that each monomer- water 
contact involves one cluster. The Hamiltonian follows 
from Eq. (p|): 



KmL 



(15) 



where N^y — {z — 2)N -|- 2 is the highest number of water 
sites that can be in contact with the polymer, and Ij is 
the number of contacts between the j-th water site and 
the monomers. The partition function of the polymer in 
solution reads as: 



Z = ^Z(C) = ^C(ric)^(nc) 



(16) 



where C are the conformations of the polymer and Z{C) 
the restricted partition function, obtained tracing over 
water variables at fixed conformation C. Due to the 
form of the Hamiltonian given in Eq. (15), Z{C) de- 



pends only on the total number of water-monomer con- 
tacts ric, and ^ (Uc) is the number of SAWs character- 
ized by the same value of ric. Z can be factorized as 
Z = ZfaZ/, where Zf, = (Z^^'")^^"' = A^^"' is the contri- 
bution of all water sites when in contact to other water, 
and Zi = C(nc)a;"'=. According to the above fac- 
torization, the specific heat is the sum of a bulk contri- 
bution Cb and of an interaction contribution C/. We in- 
troduce also the average number of water-monomer con- 
tacts, (ric) = x-g^ logZ/, which is a measure of the com- 
pactness of the polymer. To exactly evaluate the above 
quantities, Q (uc) should be obtained through an exhaus- 
tive enumeration of the SAWs. However, if we restrict 
ourselves to a square lattice, the numerical estimates re- 
ported in |l^ allow us to write 



C (n,) ~ Co^^(aoA^)"'^"=) exp(aoiV) , (17) 

where w(ric) — {Niy — nc)/2 is the number of monomer- 
monomer contacts, ao = 0.164 and Co is the number of 
SAWs of length N. Such an estimate is expected to be 
very good if N is large jl^. Hence, an analytical expres- 
sion can be found for Zj: 
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c^-r(u;„,.-Hl,^)e--^(i-:^), (18) 



where r(n, x) is the incomplete F- function and lUmax = 
{Nw — n™'")/2 (the smallest number of contacts is the 
perimeter of the globule, n™'" ~ 2\/lV7r). The results for 
Ci and (ric) are reported in Fig. ||. 

The presence of both "cold" and "warm" collapse tran- 
sitions, signalled by the drop of the number of contacts 
and by the jumps in the specific heat is strikingly evident. 
Notice that the polymer swells at low temperature even 
with K = 0: water-clusters' energy at T = is not af- 
fected by the contact with monomers, and cold unfolding 
is not a trivial effect of a swollen-ground-state dominance. 
As N grows both transitions get sharper, thus suggesting 
the existence of true phase transitions in the thermody- 
namic limit. The above phenomenology is very close to 
that of proteins in solution Moreover, swelling at 
low temperatures has been found in homopolymers like 
poly(N-isopropylacrylamide) (PNIPAM) and poly(N- 
isopropylmethacrylamide) (PNIPMAM) [^, while warm 
swelling has been detected for a number of homopolymers 
in non- aqueous solution |18|. 

Both PNIPAM and PNIPMAM present a polar peptide 
group in the side chain and hence are able to form hydro- 
gen bonds; hence, one could object that it is not appro- 
priate to propose a comparison with their phenomenol- 
ogy, since polar groups favourably interact with water 
and help the swelling of the chain. Actually, we are not 
aware of any calorimetric experiment on completely non- 
polar homopolymers in water (that would represent ideal 
tests to our predictions, provided that they swell in the 
range T = — 100 °C and the solution is so dilute that 
aggregation is negligible) . However, hydrogen bonds can 
be formed both in the coil phase, with water, and in the 
globule one (and are indeed considered a possible cause 
for the detected hysteresis in the coil-globule transition 
]lOt). Assuming a negligible overall contribution of hy- 
drogen bonds to the stability of the coil over the globule 
phase, one can admit that the cold transition in PNI- 
PAM and PNIPMAM is mainly due to the temperature 
dependence of the free-energy of hydration of the nonpo- 
lar groups. This is why we think that it is meaningful to 
compare our prediction for a nonpolar homopolymer to 
experimental results on PNIPAM and PNIPMAM. Our 
model is not intended to give a detailed description of 
their phenomenology, though: because of its dimension- 
ality, of its lack of geometric detail, and of the assump- 
tions that make it simple, it is not suited to describe the 
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"domain" cooperativity within of the coil-globule transi- 
tion, or the freezing of side chains within the globule. 

Warm swelUng of PNIPAM and PNIPMAM has not 
been detected: this is due to the fact that the tempera- 
ture of warm transition should lie above the water boiling 
temperature [plj. 



III. CONCLUDING REMARKS 

We have presented a simple model which reproduces 
the two-dimensional simulation pi and the experimental 
results [| for the trends of AF,1:,E, AS, AC upon hy- 
dration of a nonpolar solute, showing that the reduced 
number of ways of forming hydrogen bonds in the pres- 
ence of a non-bonding solute, together with a shift in 
hydrogen bonding energy, are sufficient to reproduce the 
experimental hallmarks of hydrophobic hydration. This 
is consistent with the views expressed in Refs. [^,|| and 

®- . . . 

Our description of water is focused on the clusters in 

the hydration shell, which are affected by the presence of 
the solute molecule. To keep the model analitically inte- 
grable, we disregard interaction between clusters. In this 
way, we are able to describe hydration, but cannot de- 
scribe solvent-mediated long-range interactions between 
solute molecules. For this reason, this approach is some- 
what complementary to the one recently proposed by 
Kolomeisky and Widom ||], which can be treated an- 
alytically in the one-dimensional case, and is very well 
suited to evaluate the potential of mean force. 

When applied to a homopolymer in solution, our model 
recovers the experimental swelling at both high and low 
temperatures. We reported only a calculation for a two- 
dimensional polymer. We did so for the sake of the inter- 
nal coherence of the model, because the geometric con- 
tribution to the density of states of water clusters is es- 
timated for the two-dimensional case and the entropies 
of both the water and the polymer play a key role in the 
behavior of the polymer in solution. Anyway, the fact 
that our two-dimensional water yields a correct qual- 
itative description of nonpolar hydration suggests that 
the application of the 2-d water model to the hydration 
of 3-d SAWs could be interesting; indeed, we have per- 
formed some calculations using three-dimensional SAWs 
and have verified that they also reproduce the same qual- 
itative results; of course, the temperature range of the 
compact phase varies a bit. 

The fact that our model describes both cold and warm 
swelling of a nonpolar homopolymer suggests that it 
could be effectively used in protein folding studies, de- 
scribing both cold and warm unfolding in a single frame- 
work. Even if hydrophobicity alone cannot be safely con- 
sidered as the unique responsible for protein folding and 
stability, it is very likely that any progress in its descrip- 
tion, and in its encoding in simple models, will yield signi- 
ficative improvement in the understanding of the protein 



folding puzzle. 
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FIG. 1. Cluster of water molecules. Bonds may be formed 
or broken, according to energy. The circle represents a so- 
lute (shell case) or (an)other water molecule(s) (bulk). Left: 
ordered arrangement with molecules in state "E" (one "ex- 
ternal" (e) and two "on-shell" (s) arms). When all bonds are 
formed this is the ground state. Right: in state-I molecules 
have two e arms and an internal (i) one, which cannot form 
a bond if a solute is present. 
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FIG. 2. Free energy AF (dash-dotted), energy AE 
(dashed), entropy TAS (dotted) and specific heat AC/3 
(solid line; 1/3 is for rendering purposes) changes upon hy- 
dration of a nonpolar solute. Here 2 = 4, m = 4, = 0, 
Qb = rn/2, hs/ht — 1.2, J = 0.3 (so that K = 0). Boltzmann's 
constant ks ~ 1 so that the specific heat is adimensional. En- 
ergies and temperatures are in units of hb- The shaded region 
can be compared with experiments M and simulations M. 
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FIG. 3. Nonpolar homopolymer in solution: temperature 
profiles, for different lengths A'', of the excess specific heat 
per water site, cj — Ci/Nw (left scale; dashed line: A'^ = 200; 
dotted, N — 10®), and average fraction of water-monomer 
contacts, = {nc)/{2N + 2) (right scale; dash-dotted line: 
N = 200; solid, N = 10*^). Parameters as before. 
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